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CHAPTER 1 


INTRODUCTION 

Airplane optimization has always been the goal of airplane designers. In the 
conceptual design phase, a designer’s goal could be tradeoffs between maximum 
structural integrity, minimum aerodynamic drag, or maximum stability and control, 
many times achieved separately. Bringing all of these factors into an iterative 
preliminary design procedure was time consuming, tedious, and not always accurate. 
For example, the final weight estimate would often be based upon statistical data from 
past airplanes. The new design would be classified based on gross characteristics, such 
as number of engines, wingspan, etc., to see which aiiplanes of the past most closely 
resemble<fthe new design.(Ref. 1) This procedure works well for conventional airplane 
designs, but not very well for new innovative designs. 

With the computing power of today, new methods are emerging for the 
conceptual design phase of airplanes. Using finite element methods, computational 
fluid dynamics, and other computer techniques, designers can make very accurate 
disciplinary analyses of an airplane design. These tools are computationally intensive, 
and when used repeatedly, they consume a great deal of computing time. In order to 
reduce the time required to analyze a design and still bring together all of the disciplines 
(such as structures, aerodynamics, and controls) into the analysis, simplified design 
computer analyses are linked together into one computer program. These design codes 
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are very efficient for conceptual design. One such program, called ACS YNT (Ref. 2-4), 
is being developed by the NASA Ames Research Center. 

The work in this thesis is focused on a finite element based conceptual design 
oriented structural synthesis capability (CDOSS) tailored to be linked into ACS YNT. 
CDOSS and ACS YNT exchange data. CDOSS automatically generates a mesh for a 
finite element model for airplane wings from the geometry data received from 
ACSYNT. The structural analysis of the finite element model includes stress analysis 
with calculation of deformations and buckling stresses. Optimality criteria based 
structural optimization is used in CDOSS. Weight is estimated by summing the weight 
of the individual elements corrected for manufacturing information. 

Shape variation is also very important in conceptual design. Changes in local 
elements as well as overall planform should be part of the optimization process. The 
structural analysis is done by generating a mesh for a finite element model. The 
difficulty with finite element models is that greater accuracy is obtained by using a finer 
mesh. However, a finer mesh requires more computing time. Optimization is 
performed using constraints on deformations, stresses, and buckling. The focus of this 
work is to perform these tasks efficiently and accurately for conceptual design. 

The outline of this work is as follows: in chapter 2 the finite element modeling 
and analysis is discussed. Chapter 3 focuses on the wing shape parameterization and 
configuration. In chapter 4, the behavior constraints for optimization are described. 
Chapter 5 has a discussion of the structural weight evaluation and optimization 
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algorithms. In chapter 6, the test case of a fighter type wing is discussed. Information 
concerning the use of CDOSS is in the appendix. 



CHAPTER 2 


FINITE ELEMENT WING MODELING AND ANALYSIS 

2.1 Introduction 

The modeling for wing box finite element analysis is described in this chapter. 
The construction of the wing and the elements used will be discussed. Sensitivities with 
respect to shape will be described. Finite element analysis of the wing model and 
calculation of stresses and deformations will be covered. A brief discussion of the 
structural vibration modes will also be included in this chapter. Load cases and 
interpolation of loads to the finite element mesh will be described. 

2.2 Wing Construction 

The wing box is constructed with cover skins, spars and ribs. Membrane (plane 
stress) elements and truss elements are use to model the parts of the wing box. Spars 
and ribs are modeled with two caps and a shear web between the caps. Cover skins and 
webs are modeled using membrane elements in plane stress. Truss (rod) elements are 
used for spar and rib caps. Figure 1 illustrates the wing box model. The cover skin is 
made of thin membrane elements. The spars (or ribs) are modeled with caps and a thin 
shear web. Skin membrane elements can take tension, compression, and shear stresses. 
Web membrane elements are in shear only and cannot take tensile or compressive 
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stresses. The rod (truss) elements are in tension or compression and cannot take 
bending stresses. Thus, the wing skin and spar caps resist resist transverse shear. 
Vertical spacers (rods) are placed at all nodes between the upper and lower skin. The 
spacers are necessary to keep the skins from collapsing onto each other. 



Figure 1. Wing Box Model 
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2.3 Elements 

The elements chosen are simple and allow for faster computation and closed 
form, explicit analytic sensitivity analysis. (Ref. 5) The rod elements are simple linear 
elements. The membranes are broken into constant strain triangular (CST) elements. 
All elements have constant stress throughout the element which causes stress 
discontinuities where elements meet. In parts, such as the wing skin, the stresses can be 
smoothed to get a continuous change of stress values from root to tip. For a more 
refined grid, “dummy” elements are used to support nodes which are not on actual ribs 
or spars. In figure 2, a portion of the wing with three ribs and 3 spars has 9 nodes and 
the elements are quite large. With the addition of two “dummy” ribs and two “dummy” 
spars, the same section of the wing now has 25 nodes. The length of each truss element 
is halved and the area of each membrane element is quartered. “Dummy” spars and 
“dummy” ribs are used solely to create a better, more refined mesh. They do not 
provide any structural strength. To ensure that they neither add to strength or weight, 
these elements are given very small thickness (1% of that of real elements) and we can 
add nodes without causing singularities in the stiffness matrix.(Ref. 6) The real 
elements and “dummy” elements are treated the same way in the analytic differentiation 
with respect to shape. 
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2.4 Sensitivities 

Structural behavior sensitivities with respect to shape design variables is 
important for shape optimization. Sensitivities are used to calculate gradients and to 
construct approximate constraint functions and objective functions. Since we used 
simple finite elements in modeling the wing, the mass and stiffness matrices are explicit 
and algebraic in nature. This allows us to differentiate the basic structural equations 
explicitly with respect to any design variable. We obtain sensitivities in a closed, 
explicit analytic form without numerical integration. (Ref. 7) 

2.5 Deformations 

For a static structural system, the governing equation is given by 

[K]{u} = {F} 

where [K] is the banded global stiffness matrix, {u} is the global displacements vector, 
and {F} is the nodal loads vector. We solve for the displacement vector using a 
decomposition technique and skyline solver of Ref. 8. 

2.6 Stresses 

Individual finite element stresses are calculated using the global displacements 
found previously. For rod elements, we use Hook’s stress/strain law to find the axial 
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stress in the element. (Ref. 9) The change in the length of the rod is found from 
subtracting the displacements of the two ends. Strain is the change in length divided by 
the original length of the rod. Then using Hook’s law, a=Ee, we find element stress. 
The local strains can be transformed into global strains in order to obtain global stresses. 
For the CST membrane elements, the procedure is similar. The stress is no longer a 
scalar. It is a vector of a„, a yy , and a xy and E is a matrix of constitutive properties. 

2.7 Vibration Modes 

The vibration modes are important when analyzing a dynamic structural system. 
In the case of wing structures, they are required for flutter and dynamic response 
analysis. To find the natural frequencies and the vibration modes, we use an eigenvalue 
function method. The governing equation for undamped simple harmonic motion is 

[K - orM] {<|>} = {0} 

K is the global stiffness matrix, co is the natural frequency, M is the global mass matrix, 
and <{> is the mode shape. For a non-trivial solution, the determinant of [K - co 2 M] must 
be equal to zero. The eigenvalues of the system are equal to © 2 and corresponding 
eigenvectors are mode shapes. 
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2.8 Loads 

There is a large number of loads that can be applied to a wing. These loads 
come from many different sources. Some possibilities are lumped mass from external 
stores or tanks, uniform pressure loads, distributed loads, loads from control surface, or 
aerodynamic loads. These loads must be applied carefully to the planform of the finite 
element model. For a finite element model, all loads must be applied at the nodes. 
Figure 3 illustrates how a uniformly distributed load is transformed into concentrated 
loads at the points where there is structural support. Nodes where “dummy” spars and 
“dummy” ribs intersect are called “floating” nodes. Loads cannot be applied at 
“floating” nodes where no real structure exists to support the load. We must use 
interpolation to “move” loads to nodal points while keeping the load distribution and 
force resultants as close to the original as possible. 



Distributed Load (Ih/in) 



Interpolated Load (lbs) 

(Equivalent concentrated loads at nodal points) 


Figure 3. Load Transformation 
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2.9 Aerodynamic Load Transformation 

Aerodynamic analysis of a wing utilizes a grid of points similar to the grid used 
in the structural analysis. Usually the points of these two grids will not match one 
another. In order to use the aerodynamic load distribution in the structural analysis, the 
loads must be transformed to the structural grid. One method is to use energy principles 
and the concept of work. Deformations can be represented by a series of polynomial 
functions or by shape functions. A description of both types of function transformation 
follow. A comparison of results using the Chebychev functions and the shape functions 
will be presented in chapter 6. 

2.9.1 Polynomial Functions 

Displacements are represented by a series of polynomial functions. 

w(x,y) = f 1 (x,y)qi+f 2 (x,y)q 2 + . . . +f„(x,y)q n 
The series is evaluated at all structural points and aerodynamic points. All the grid 
points are transformed to a unit box to help prevent ill-conditioning in matrices. The 
two matrices, [Ws] and [Wa], can be created where each row of [Ws] consists of the 
polynomials, f„, evaluated at a structural point and each row of [Wa] consists of the 
polynomials evaluated at an aerodynamic point. Requiring that the work done by each 
system is equal results in the following matrix equation. 
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{ Fs } T [Ws] {q } = { Fa} T [Wa] {q } 

{ Fs } and {Fa} are vectors of the forces at the structural points and aerodynamic points 
respectively. The equation above must hold for any vector {q}. The following 
derivation results in a relationship to find {Fs}, the equivalent forces on the structural 
grid. 

{Fs} T [Ws]= {Fa} T [Wa] 

[Ws] T {Fs} = [Wa] T {Fa} 

[Ws][Ws] t {Fs} = [Ws][Wa] T {Fa} 

Let 

[A] = [Ws][Ws] T 

[B] = [Ws][Wa] T 

where [A] is now a square, symmetric matrix. 

[A]{Fs} = [B]{Fa} 

{Fs} = [A] 1 [B]{Fa} 

Using the above equation, the equivalent force on the structural points can be found. 

The polynomial functions chosen for this work were Chebychev polynomials. 
The Chebcychev polynomials are formed in the following manner. 

Ti(x) = 1.0 
T 2 (x) = x 
T 3 (x) = 2x 2 -1 
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T„(x) = 2xT n .|(x)-T„. 2 (x) 

The polynomial functions are evaluated for the x and y coordinates for both the 
structural grid points and the aerodynamic grid points. Then the matrices [Ws] and 
[Wa] can be created, where a typical shape function f r (x,y) is given by 

fr(x,y) = Tj(x)*Tj(y) 

Chebychev polynomials use high order terms, but do not have the same ill-conditioning 
problems of simpler polynomials, such as f n (x) = x n . 

There are a few limitations when using a polynomial series to represent the 
displacements of the entire wing. First, since the grid is transformed into a unit box, the 
points must be somewhat uniform on the grid. For instance, if there are 10 points on the 
root chord, then there should be 10 points on each successive chord. Second, the 
number of terms in the polynomial series must be equal or greater then the number of 
structural points. [A] = [Ws][Ws] T is singular if this condition is not met. 
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2.9.2 Shape Functions 

Shape functions can be used to represent deformations in functional form. 
(Ref. 10) Deformations of a structural set of points can be interpolated to an 
aerodynamic set of points. The grid of structural points is first transformed into a 
rectangular grid as illustrated in figure 4. Each square cell is an isoparametric eight- 
noded quadrilateral. (Ref. 1 1) Only nodes that lie on a real spar or rib are of interest. 
We use the following quadratic interpolation functions on each cell: 

JV,=-7(l-5)(l-0(l + s + 0 
4 

(V 2 = - t(1 + ■*)(! - 0(1 -s + t) 

4 

^ 3 =-7(l + s)(l + r)(l-5-0 
4 

W 4 =-^-(1- j )(1 + 0(1 + 5-0 
Nj =“(l-J 2 )(l-0 

N 6 =^(1 + *)( 1 -' 2 ) 
w 7 =^(i-* 2 )(i + 0 

W g =^(l- 5 )(l-f 2 ) 

The coordinates s and t are normalized coordinates for each cell as shown in figure 4. 
By locating the cell in which an aerodynamic point lies, the displacement of the 
aerodynamic point can be found when performing the following summation over the 


structural points of the cell. 
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*.=Zz,W. 

1*1 

The transformation matrix, D, can be found by cycling through all aerodynamic points. 

Z a =D9 

where Za and 9 are vectors of displacements for the aerodynamic points and the 
structural points respectively. By the following derivation, the proper transformation for 
forces is found. Begin by requiring that the work done by forces in the aerodynamic 
grid, F a , be equal to the work done by the forces on the structural grid, F s . 

FJ 9 = F T Z 

s a a 

Fj9 = F/D9 

F T = F t D 

* j a 

F = D T F 

A s a 


Using this relationship to find equivalent forces on the structural grid results in forces 
that preserve the conservation of forces and moments. 
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Figure 4. Grid Transformation 





CHAPTER 3 


WING SHAPE PARAMETERIZATION 


3.1 Introduction 

This chapter describes how the model is constructed. It is important to model 
the shape of the wing box as accurately as possible. The wing shape geometry is defined 
by ACS YNT. The automatic generation of the finite element mesh will be discussed. 
The internal geometric layout of spars and ribs is also covered in this chapter. 

3.2 Wing Segmentation 

The wing is segmented into trapezoidal sections to help maintain a more accurate 
definition of the shape of the wing box. Each section has a constant sweep angle. A 
change of sweep would lead to a new section. Figure 5 illustrates the segmentation of 
the planform of a Concorde wing. This allows the planform of each section to be a 
simple trapezoid which can then be easily divided into triangular membrane elements. 
Segmenting the wing in this manner allows for great flexibility in new innovative 
planform designs. Partially swept forward wings can be tested and models for wings, 
such as the Concorde wing, can be easily and accurately constructed. 
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Figure 5. Wing Segmentation 





3.3 ACSYNT Geometry Definition 


The geometry of the wing is obtained from ACSYNT. Following the same 
segmentation scheme as previously described, ACSYNT provides airfoil cross-sectional 
data for each segment. Beginning at the root chord and continuing to the tip chord, 
cross-sectional data consisting of the airfoil shape is given for each segment. The airfoil 
shape is defined by xyz coordinates of a point on the airfoil. These points are linearly 
connected from root to tip of a segment to form the shape of the airfoil. The same 
number of points are given for each cross-section and the listing of the points always 
begins at the trailing edge. This creates a one-to-one correspondence of points from one 
cross-section to the next. These points are connected linearly to form the upper and 
lower surface of the wing. Figure 6 shows a typical segemented wing geometry. The 
points on the airfoil that are data points from ACSYNT are shown. Also illustrated in 
figure 6 is the manner in which CDOSS connects the points to form the wing surfaces. 
Since the leading edge and trailing edge control surfaces do not provide significant 
structural strength, we exclude them and only analyze the wing box. The control 
surface are defined as a percent of the chord. The leading edge and trailing edge spars 
are placed at the leading edge of a trailing edge flap and at the trailing edge of a leading 
edge flap, respectively, to connect the upper and lower surface. The wing box is 
constructed and ready for mesh generation. 
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Figure 6. Wing Geometry Using ACSYNT Data 
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3.4 Automated Finite Element Model Generation 

The capability for automated finite element model generation is essential for 
efficient conceptual design analysis. (Ref. 12) Creating finite element input files by 
hand is both time-consuming and tedious. CDOSS has the capability to take the 
geometry given by ACSYNT and create a finite element model. Automatically, all of 
the nodes, membrane elements and truss elements are created automatically using the 
following guidelines. Nodal points exist on the upper and lower surface anywhere spars 
and ribs intersect (including “dummy” spars and ribs). These intersections also define 
quadrilateral cells on the skins, each with four nodes. These quadrilateral cells are 
divided into two triangular membrane elements. Webs for all spars and ribs are divided 
into quadrilateral cells between upper and lower nodes along the spar and rib lines. 
These quadrilateral cells are also divided into two triangular elements. Truss elements 
are placed between nodes along spar and rib lines on the upper and lower surface for the 
real structural elements only. “Dummy” spars and ribs do not have caps. All elements 
are also assigned thickness or area while generating the elements. Figure 7 illustrates 
how all of these elements are used to model a wing using the automatic mesh generator. 
These guidelines make it possible to generate the entire wing model where nodes are 
defined by geometric coordinates, elements are defined by nodes, and design variables 


are assigned appropriately. 



Caps 


Vertical spacers 


Spar or Rib 


Figure 7. Mesh Generator Example 
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3.5 Layout of Spars and Ribs 

For the conceptual design phase, flexibility in the design is important. One 
aspect of the design is the layout of spars and ribs in the wing box. CDOSS is capable 
of two different spar configurations. The first option is ‘fanned’ spars. Spars are spaced 
evenly at each cross-section and extend from the root to the tip. The spars are 
continuous with bends at section breaks. The second option is parallel spars. In this 
configuration, the spars are parallel to the trailing edge of the wing box. The spars are 
evenly spaced at the root chord and run parallel to the trailing edge of each section. The 
spars are continuous to the tip unless they intersect the leading edge, where the spar 
would then end. At each section break, the spar continues but runs parallel to the 
trailing edge of the next section. Figure 8 illustrates both of these configurations. In 
both spar configurations, the number of spars beginning at the root is a design variable. 
“Dummy” spars are added between real spars to refine the mesh and follow the same 
configuration guidelines as real spars. 

The layout of ribs depends on the choice of spar configuration. CDOSS is 
presently limited to ribs parallel to the root rib. The number of ribs is a design variable. 
Each section can have a different number of ribs. For fanned spars, the ribs are spaced 
evenly between section breaks and “dummy” ribs are added between real ribs to refine 
the mesh. For parallel spars, ribs are placed at points in which a spar intersects the 
leading edge. This is necessary to ensure that quadrilateral cells are created for finite 
element model generation. If additional ribs are desired, they are placed between 
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CHAPTER 4 


BEHAVIOR CONSTRAINTS 


4.1 Introduction 

Behavior constraints are the key to optimization. The constraints represent 
design limitations and failure modes. Deformation constraints may represent stiffness 
requirements. Stress constraints protect against failure in yield or fatigue. We use a 
simple comparison of the stress in an element to the maximum yield stress, taking into 
account the margin of safety. Buckling constraints are very critical for wings. As the 
wing bends up in flight, the upper surface experiences large compressive stresses. 
Behavior constraints are discussed in detail in this chapter. 

4.2 Deformation Constraints 
4.2.1 Boundary Conditions 

Wing deformations are constrained as follows. The root chord is fixed to a 
reference surface. All degrees of freedom are not allowed to move at the root chord. 
The rotational degrees of freedom are fixed throughout the wing. This is because of the 
use of CST and truss elements which cannot resist local rotation at their nodes. 
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4.2.2 Behavior Constraints 

Constraints on the amount of deformation at the nodes or deformations of the 
elements in relation to each other could also be constrained. 

4.3 Stress Constraints 

Each element is associated with a behavior constraint on stress. The material 
has a given yield strength that each element cannot violate. The truss elements have a 
single, scalar axial stress. The axial stress multiplied by the margin of safety cannot be 
greater than the yield strength. We have assumed that compressive yield strength is the 
same as tensile yield strength. Membrane elements have two in-plane stresses and a 
shear stress. We use an equivalent Von Mises stress for comparison with the maximum 
yield stress. 

Equivalent Stress = sqrtfcti 2 - G1O2 + C2 2 ) 

where <7i and Gi are principle stresses. The equivalent stress multiplied by the margin 
of safety cannot be greater than the yield strength. 
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4.4 Buckling Constraints 

Panel buckling constraints are very important in the optimization routine of 
airplane wings. In most cases, the buckling constraint is the determining factor for 
changes in the design variable for the wing skin. In order to perform a buckling 
analysis, the panels must first be identified and panel stresses calculated. Then, 
buckling constraints can be evaluated. 

4.4.1 Panel Identification 

Panels must be defined by the real structure where there is rib or spar support. 
Intersection of real spars and ribs serve as supports for the panels on the wing skins and 
in the spar and rib webs. By defining the panels this way, each panel has many 
triangular membrane elements within it. In figure 9, the section of the wing skin shown 
has three real spars, two “dummy” spars, three real ribs, and three “dummy” ribs. This 
section has four buckling panels. Each panel has eight triangular membrane elements. 
The dummy elements do not provide support for the buckling panels. Each panel is 
identified by the nodes at the four vertices and the membranes that make up the panel. 
The average of the thickness of the membranes in a panel is used as the panel thickness. 
Using the coordinates of the four vertices, the edge lengths are calculated. The panel 
dimensions are the average of the two lengths in each direction. Thus, each trapezoid 
panel is approximated as a rectangular panel. Membrane normal and shear stresses are 
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averaged to obtain the panel normal and shear stresses. These dimensions and stresses 
are found for each panel, so that simplified buckling analysis can be performed on each 
panel. 



Figure 9. Panel Identification 
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4.4.2 Buckling Analysis 

A simple conservative approach was adopted for the buckling analysis. We 
assume each panel to have simply supported edges. Actually, there is always more 
rotational stiffness along panel edges. This simplified buckling analysis used ratios of 
the panel stresses to the critical buckling stresses. Stresses are rotated to local stresses 
since the rectangular panel is oriented in local coordinates. Figure 10 illustrates the 
approximation of the trapezoidal panel as a rectangular panel. The local stresses used to 
perform the buckling analysis are also shown. Critical buckling stresses are functions of 
material properties and panel dimensions. (Ref. 13-15) A simple interaction equation 
was used for the constraint. Since the normal stresses and shear stresses all affect one 
another, we cannot perform the buckling analysis for each stress individually. The 
buckling interaction equation is as follows (Ref. 17): 

Nz Ny (Nxy V 

g ~ N*, Ny~ 

Critical buckling stresses are negative. Compressive panel stresses are also negative. If 
g, the constraint value, is negative, then the panel will buckle. This procedure is carried 
out for each panel to check for buckling. Another approach to buckling analysis is a 
global approach where constraints are based on a linear eigenvalue analysis of the entire 
structure. (Ref. 18) The advantage of this approach is that it allows interaction of the 
individual panels. The disadvantage is additional computation time is needed to solve 
an eigenvalue problem. The simplified local panel buckling analysis used in CDOSS is 
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much faster. Also, this approach is easily applied to composites by simply modifying 
the equations for panel stresses and critical stresses. (Ref. 19) 



Figure 10. Panel Approximation 




CHAPTER 5 


STRUCTURAL WEIGHT EVALUATION AND OPTIMIZATION 

5.1 Introduction 

It is very difficult to obtain accurate wing weight for the conceptual design 
phase. Details, which are not a part of conceptual design, account for a considerable 
portion of the weight. Such details are fittings, joints, attachments, and access doors, 
for example. There are different methods of estimating structural weight which will be 
discussed in the chapter. We will also cover the optimization algorithms used. CDOSS 
uses a simple optimality criterion, although there are more rigorous nonlinear 
optimization techniques which require more computation time. (Ref. 20-25) In this 
chapter, we will discuss the technique used in CDOSS. 

5.2 Weight Evaluation 

Several different methods of weight estimation have been used in conceptual 
design. One method is statistical weight estimation. (Ref. I) This method is dependent 
on a database of existing airplane designs. Statistical weight estimation does not give 
accurate results for non-conventional airplane designs and does not take into 
consideration more detailed design variations such as chord extensions or varied wing 
thickness distribution. (Ref. 26) Another method is a theoretical estimate of bending 
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weight. (Ref. 26) Many constraints are accounted for using this method, such as 
bending stress, buckling and local pressure loads. A more structurally efficient design 
can be achieved. The disadvantage of bending weight estimates is that the results are 
less accurate as particular considerations, such as buckling, dominate the sizing of the 
structure. (Ref. 26) 

CDOSS uses a method that is more accurate than statistical estimates, simpler 
than bending weight, and takes advantage of the data calculated in the finite element 
analysis and optimization. All of the major structural elements are identified for the 
finite element model. All of the elements of the finite element model are sized within 
the optimization routine. To obtain the structural weight of an airplane wing, we simply 
calculate the volume of each element, sum the volumes, and multiply by the density of 
the material. This method results in fairly accurate structural weight estimates for 
conceptual designs using little computing time. Any approximations in modeling the 
structure will change the final weight estimate. Statistical, empirical correction factors 
used to convert the idealized weights of skins, caps and webs to real weights taking 
manufacturing techniques into account. The structural weight estimating technique used 
in CDOSS is an efficient method for preliminary design. 



5.3 Structural Optimization 


34 


Automated multidisciplinary optimization is necessary to numerically handle 
many competing design constraints. (Ref. 27) While rigorous mathematical 
programming techniques are established for structural synthesis (Ref. 24), the method 
used in this work is an optimality criterion method using fully-stressed design and 
buckling constraints for structural optimization. This is done mainly for the sake of 
simplicity at this stage. (Ref. 29) 

5.3.1 Element Sizing 

This is the approach used in CDOSS. A fully-stressed design is based on stress 
constraints. Each individual element is sized based on a ratio of element stress to 
maximum yield stress. (Ref. 29) The thickness or area of the element is increased or 
decreased based on the ratio of the element stress multiplied by a safety factor and of the 
maximum yield. The new size is the ratio multiplied by the present size. Since strength 
is a linear function of thickness stress (assuming internal forces do not change much 
from iteration to iteration), this method is effective. Analytic sizing for buckling is 
much more difficult. Buckling stresses are cubic functions of thickness. Within the 
constraint equation, the in-plane normal stresses are not of the same order as the shear 
stress. (Ref. 21,24) For these reasons, the sizing for buckling is based on a pre-selected 
value of the buckling constraint. The decision of whether to increase or decrease the 
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thickness of the panel is dependent on the value of the constraint function. If the 
constraint value is negative, the panel is buckling and the thickness should be increased. 
If the constraint value is positive, then the thickness should be decreased. The sizing for 
buckling iterates until all panels have a constraint value between zero and one tenth. 

5.3.2 Optimization Routine 

The optimization routine calculates the new thickness for stress and buckling 
constraints separately, then compares them. The new thickness or area does have 
maximum and minimum limits that are determined from other constraints, such as 
manufacturing constraints. The new thickness from the stress and buckling constraints 
are compared to the minimum limit. The largest of these three values is compared to 
the maximum limit. The smallest of these two values is the new thickness assigned to 
the element. This optimization routine will iterate until the solution converges. 
Following each iteration, stresses and weight are calculated based on the new values for 
thickness. The convergence check first requires that no stress constraint is violated by 
more than 5 %. The second requirement for convergence is that the percent change in 
weight between the present weight and the weight of the previous iteration is not more 
than a selected tolerance criteria. The method of optimization is effective in ensuring 
that both stress and buckling constraints are not violated. Since the calculations are 
simple, this optimization routine is also very efficient. 



CHAPTER 6 


TEST CASE: FIGHTER TYPE WING 


6.1 Introduction 

The test case used for this work was a typical fighter all aluminum wing. This 
chapter discusses the modeling of the wing and the load cases used. An analysis of the 
load results will also be covered in this chapter. Optimization and weight results will 
also be discussed. All measurements are in English units (inches, pounds, etc.). 

6.2 CDOSS Wing Model 

The CDOSS fighter wing model very closely resembles an actual wing. The 
CDOSS model is a model generated from ACSYNT cross-section data using the design 
options available in CDOSS. The table 1 shows a comparison of major features of the 
CDOSS model and an industry quality finite element model. 
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Table 1. Comparison of CDOSS and Industry Models 


Feature 

CDOSS Model 

Industry Model 

Number of Spars 

12 

1 1 

Number of Ribs 

5 

5 

Length of Root Chord 

85 

80 

Length of Tip Chord 

20 

10 

Half Span 

140 

130 


The spar and rib layout match nearly exactly with spars parallel to the trailing edge and 
ribs parallel to the root chord. The dimensions are only slightly different. The CDOSS 
model is a little larger than the actual wing box. Figures 1 1 and 12 show the finite 
element meshes for the CDOSS model and the industry model. No “dummy” spars 
were used for the CDOSS model of the wing. There are five “dummy” ribs, represented 
as dotted lines in figure 11. There are a few other minor differences in the finite 
element meshes, such as more structure and more detailed mesh on the industry model 
where pylons or external stores would be supported. Also, the industry mesh used 
quadrilateral elements rather than triangular elements. CDOSS automatically generated 
a very accurate model of the fighter wing box. 
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6.3 Load Cases 

Two load cases were analyzed for the fighter wing. The first load case was a 
simple 1000 pound point load applied to the rear spar at the tip. The second load case 
was a 9 'g' pullup at .9 mach, 10,000 feet altitude, 21,289 pounds gross weight, and tip 
missiles. Load data was given as concentrated loads at the nodes of the industry mesh. 
These loads were interpolated to the nodes of the CDOSS mesh. Since there are no 
“dummy” spars, there are no “floating” nodes. Thus, every node on the mesh is 
supported by real structure and can carry loads. Since the meshes are so closely 
matched, the interpolation was very simple. The tip load is a simple load case that is 
not representative of realistic loads on a fighter. It is possible to use this load to predict 
the displacements and reaction of the model and to assess convergence characteristics 
and effects of various constraints. Also, empirical data was available for the industry 
model and from an Automated Structural Optimization System(ASTROS). (ASTROS 
is being developed by the U.S. Air Force.) The ASTROS model has the same mesh as 
the industry model. The 9 ‘g’ load is a typical maneuver of a fighter airplane. This load 
case models one of the critical maneuvers which the airplane is designed to perform. 
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6.3.1 Comparison of Transformation Methods 

In the case of the aerodynamics loads, forces must be transformed to the 
structural grid. The two methods described in chapter 2 used functions to represent the 
deformations. The shape function method was applied to individual cells. The 
Chebychev polynomial method was applied to the entire wing. Since the Chebychev 
method requires that the structural grid be uniform, a different model was used for this 
comparison. The model is illustrated in figure 13. There are 5 real spars and 6 
‘dummy’ spars between the leading and trailing edge spars. Spars are fanned and run 
from root to tip. There are 3 real ribs and 4 ‘dummy’ ribs between the root and the tip. 
The most significant difference between this model and the fighter model used in this 
chapter is that the buckling panels are much larger in this comparison model. The 
comparison was performed using the 9 ‘g’ load case. In figure 14, the change of the 
weight per iteration is shown. Figure 15 and 16 show the wing deformation and wing 
skin thickness along a center spar line after optimization. This comparison shows that 
the two methods produce very similar results. 
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Figure 14. Comparison of Iterations vs. Weight 






% Span 


Figure 16. Comparison of Wing Skin Thickness 






6.4 Analysis of Results 
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The analysis of the CDOSS model included deformation analysis and stress 
analysis. The deformations of each nodal point is given as well as the stress in each 
element. Deformations can be used to check the flexibility of the structure. The 
stresses can be compared to the yield stress for each element. This information 
indicates whether or not the design is feasible. The deformations and stresses to be 
discussed are results following optimization. 

6.4.1 Deformations 

The deformations resulting from both the 1000 pound tip load at the rear spar 
and the 9 ‘g’ load were analyzed. A 1000 pound tip load caused a tip deflection of 
approximately 1.1 inches. This is fairly close to the 1.3 inch tip deflection of the 
industry model and the 1.01 inch tip deflection of the ASTROS model. Figure 17 
shows a comparison of the deformations from these three models. These displacements 
are along the rear spar. The CDOSS model is less flexible. There is less bending in the 
spar compared to the other models. This could be caused by the differences in spar web 
thickness and spar cap areas. The industry and ASTROS models have spar webs and 
caps with different thickness and area, respectively. CDOSS (to simplify for conceptual 
design) averages these values and gives all spar webs and caps the same averaged 
thickness or area. This will effect the reaction of the structure slightly. The 9 ‘g’ load 
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case also resulted in feasible results. See figure 18. The 8 to 9 inch deflections at the 
tip are quite reasonable and are comparable to the results of the ASTROS model. The 
analysis of the deformations of the CDOSS model showed that the model is successful 
in representing the fighter wing. 



Figure 17. Deformation Under 1000 lb. Tip Load 
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Figure 18. Deformation Under 9 ‘g’ Load 
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6.4.2 Stresses 

The analysis of stresses for both the 1000 pound tip load and the 9 ‘g’ load were 
performed. Figures 19 and 20 show that neither load case was close to the yield stress. 
The stresses plotted are taken at the mid-chord line along the span. For the tip load, the 
absolute values of normal stresses were below 2000 psi and the absolute values of shear 
stresses were below 4000 psi. These low values show that the tip load is not a critical 
load when designing this wing. The plots are fairly symmetric. When either the upper 
or lower surface is in tension at a particular point along the span, the other surface is in 
compression. The values of tension and compression are comparable and tend to 
increase from root to tip. The 9 ‘g’ load case resulted in much higher values of stress 
than the tip load. The stress was as high as 35,600 psi at one point. The value is far 
from the yield stress of aluminum (60,000 psi) and several 1000 psi less than the margin 
of safety yield stress (43,000 psi for a 1.5 margin of safety). As expected, when one 
surface is in tension, the other is in compression as the same point along the span. The 
values are not nearly as symmetric as they were for the tip load. The values tend to 
decrease from root to tip. By analyzing the stress distribution, it is possible to identify 
which sections of the wing may need extra material. 



sigma y 
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Figure 19. Stresses for 1000 lb. Tip Load 











6.5 Optimization Results 
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The optimization routine optimizes the wingskin thickness under deformation, 
stress, and buckling constraints. The results of the optimization were wingskins thinner 
than the actual fighter wingskin. Figure 21 shows a comparison of the thickness of the 
wingskin at the mid<hord line along the span. The 1000 pound tip load case produced 
extremely thin skins. This result is not surprising because this load case does not 
represent a realistic loading. Since the loading did not produce a great deal of stress in 
the structure, as was seen in the previous section. The buckling constraint was the 
determining factor for the thickness. Without the buckling constraint, the skin would 
have reduced to the minimum gage allowed without violating any stress constraints. 

This is also true for the 9 ‘g’ load case. (Minimum gage in these test cases was .005 
inches.) The wingskin thickness more closely resembled the actual wing for the 9 ‘g’ 
load because this a more critical load case. The lower surface is mostly in tension so the 
buckling constraint is not as dominant in the skin sizing. Figure 21 shows that the real 
skin thickness is much greater than thickness of the 9 g loading for the real wing. Many 
other considerations such as numerous other aircraft maneuvers, flutter, fatigue, internal 
pressure, external loads (engines, stores), etc. must be taken into account. Also, 
manufacturing constraints probably have a dominate role in the choice of the minimum 
gage This sizes parts of the wing that have little stress. Thus, the comparisons 
presented here are only for the sake of assessing global characteristics and 
computational efficiency of CDOSS. 



Upper Surface Thickness 
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Figure 21. Wing Thickness Comparison 




6.6 Weight Results 
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CDOSS calculates the weight of the wing after each analysis. When optimizing, 
the weight is calculated after each iteration. Figure 22 shows the change in the weight 
of the Fighter model as the analysis iterates until it converges. The Final weight resulting 
from the tip load is much lower than the final weight from the 9 'g' load case. This is 
the expected outcome since the tip load is not a realistic load for a Fighter. The final 
optimized weight of the Fighter model under a 9 'g' load is 416 pounds. This is a very 
accurate estimate of the structural weight of the real Fighter wing compared to 41 1.36 
pounds calculated from the ASTROS model. It should be noted that no data was 
available as to exactly what allowables were used in the industry model. Thus, several 
yield stresses were tried. Figure 23 shows that the sensitivity is fairly high at the lower 
values of yield stress. Weight is much less sensitive at higher yield stresses. The yield 
stress of aluminum alloys is approximately 60,000 psi. A material with the same 
density, but lower yield strength would result in a large increase in weight. The weight 
results for the fighter test case are excellent estimates of the real wing, thus other 
designs tested using CDOSS will also result in good estimates of structural wing weight. 
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Figure 22. Iterations vs. Weight 



Figure 23. Yield Stress vs. Weight 





6.7 Computation Time 
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CDOSS performs the structural analysis and optimization very quickly. The test 
cases for the fighter wing required only 33.6 seconds of cpu time to create a mesh, 
identify panels, and perform 5 iterations of stress analysis and optimization. (These test 
cases were run on an HP Apollo 700 UNIX based workstation.) Table 2 shows the 
breakdown of computation time for various steps of the analysis. Creating the finite 
element mesh and identifying the panels require very little time. Creating the stiffness 
matrix requires almost 1.5 seconds. The stress analysis requires 7.3 seconds. This 
includes calculation of deformations, matrix inversion using a LU decomposition 
technique, stress calculations, and a great deal of matrix manipulation. To perform one 
iteration of stress analysis and optimization requires almost the same amount of time as 
a stress analysis alone. Thus the optimization portion is extremely fast. (These time 
estimates are dependent on the size of the mesh.) Being able to perform the structural 
analysis in such a short amount of time allows designers to test many designs and 
modification quickly and efficiently. 



55 


Table 2. Breakdown of Computation Time 


Step in Analysis CPU time fsec.i 

Get data from ACSYNT & Create Mesh . 1 

Panel Identification . 1 

Compute Stiffness Matrix 1 .4 

Perform Stress Analysis 7.3 

Perform One Optimization Sequence and a Stress Analysis 7.4 



CHAPTER 7 


CONCLUSION 

The conceptual design oriented structural synthesis capability discussed in this 
thesis is tailored towards conceptual design. The analysis is accomplished efficiently 
and the results show good accuracy. The automated mesh generator quickly creates a 
finite element model for analysis using membranes elements, truss elements, and 
"dummy" elements. The simple elements allow for easy calculation, manipulation, and 
sensitivity analysis. The wing modeling is done using airfoil cross-sectional data, wing 
segmentation, and options in spar configuration. The behavior constraints used were 
effective for optimization. There are constraints on deformations, stresses, and 
buckling. Buckling is very important in the sizing of the wing skins. Using a simple 
interaction equation and averaging the stresses over panels led to simple buckling 
constraint values that could be used for optimization. 

Optimality criterion using full-stressed design and buckling constraints was 
found to be very efficient and produced adequate results for conceptual design. The 
weight evaluation is based on the finite element model corrected for manufacturing 
detail. The test case of a fighter wing showed that CDOSS modeled the wing very well 
and performed the structural analysis and optimization effectively. The results were 
compared to other models and to the actual wing. The deformations and stressed were 
analyzed. The optimization of the wing skins led to reliable weight results. 
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Computation times are extremely fast. The test case showed that the structural analysis 
and optimization of a wing in conceptual design can be cost effective and reliable. 

This work can also be extended to modeling and analyzing of an entire airplane. 
Other possible extensions for this work are in the areas of composite materials and 


added aeroelastic constraints. 
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APPENDIX A 

CDOSS PROGRAM INFORMATION 


A. I Introduction 

This appendix contains user information for CDOSS, the code developed and 
used in this thesis. The subroutines are briefly described and grouped by function. The 
program structure is outlined and main program variable are described. All programs 
and subroutines are written in FORTRAN. 

A.2 Program Subroutines 

i 

Intfem is the main program. Geometry of the wing is defined here. 

Mesh Generation: 

Interpol - interpolates points between ACSYNT provided points 
Renum - renumbers the points for convenience 
Element - creates elements for finite element model 
Shape - transforms aerodynamic loads using shape functions 

Cheb - transforms aerodynamic loads using chebychev functions 

Fload - calculates equivalent point loads for uniformly distributed loads or 
tip loads 

Panel Generation: 

Panelid - identifies the panel nodes and membranes 
Panelch - finds the dimensions and sweep of each panel 
Panelst - calculates panel average thickness, stress, and nodal forces 
FI fern is the subroutine that performs the analysis and optimization. 

The following subroutines are called within FI fern. 

Getdata is a subroutine that reads finite element data (nodes and elements). 
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Stiffness Matrix: 

Ske - computes elemental stiffness matrix for truss elements 

Kmerge - adds truss element stiffness matrix to global stiffness matrix 
Skm - computes elemental stiffness matrix for membrane elements 
Lmerge - adds membrane element stiffness matrix to global stiffness matrix 
Skyfac/Skysol solve linear equations using the method of reference 4. 

Elemental Stress: 

Stre - computes truss element stress 

Strm - computes membrane element stresses 
Buckling: 

Transform - transforms global stresses to local coordinates 
Critical - calculates buckling constraint values for all constrained elements 

Critopt - calculates buckling stress and constraints for wing skin only 

Geom calculates the geometric characteristics of the elements. 

Final outputs the results of the optimized model. 

There are also small subroutines that perform basic matrix manipulations. 

A.3 Program Structure 

A brief description of the structure of the program and the major variable within 
each section will be described. 

Wing Geometry and Mesh Generation: 

In this section, data is gathered from ACSYNT and from the user for wing geometry. 
The finite element mesh is then generated. 

Gather Data: 

herm(comp#, section#, point#,direction) - cross sectional data from ACSYNT 
sparcon - choice of spar configuration 
nspar - # of spars at the root 
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nrib(section#) - # of ribs in the section 

bshape - wing box shape 

lp,tp - % chord for slats and flaps 

sthick.rthick - spar and rib web thickness 

areal ,area2 - spar and rib cap area 

Create Elements: 

nnx - # of nodes in the x direction 

nny - # of nodes in the y direction 

upsurf - # of membrane elements on the upper surface 

xc,yc,zc(node#) - coordinates of nodes 

nnodes - total # of nodes 

mem(vertice#, membrane#) - membrane node data 
node(vertice#, truss#) - truss node dataa 
thfmembrane#) - membrane thickness 
ax(element#) - truss element area 
force(direction,node#) - applied nodal force 
ifix(direction,node#) - boundary conditions 
dum(dummy node#) - node number of a dummy node 

Panel Identification: 

In this section, panels are identified and panel characteristics are calculated. 
paneln(panel#,vertice#) - nodes at panel vertices 
panelm(panel3, membrane#) - membranes in each panel 
npanel(3) - # of panels in wingskin, rib webs, or spar webs 
panelth(panel#) - average thickness of panel elements 
pstress(panel#, direction) - average stresses in panel 
nforce(panel#,node#, direction) - nodal forces 
sweep(panel#) - sweep of panel 
apan,bpan(panel#) - dimensions of panel 
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Material Properties: 

In this section, the user inputs the material propertied needed for analysis, 
ex - elastic modulus for truss elements 
rhol - density of truss elements 
amin.amax - minimum and maximum allowable area 
em - elastic modulus for membrane elements 
rho2 - density of membrane elements 
pr - Poisson's Ratio 

tmin.tmax - minimum and maximum allowable thickness 

sf - safety factor 

yield - yield strength 

tol - convergence criteria tolerance 

Stress Analysis: 

Here the stiffness matrix is formed and the stresses calculated, 
skele - truss elemental stiffness matrix 
skmem - membrane elemental stiffness matrix 
stifk - global stiffness matrix 
disp(direction,node#) - nodal displacements 
str(truss#) - truss element stress 

smem(direction,membrane#) - membrane element stresses 
force2(node#, direction) - local forces 

Optimization: 

In this section, the optimization if performed, weight is calculated and convergence is 

checked. 

iter - iteration number 

cstrl - buckling constraint value 
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equiv - equivalent stress 

tratio - ratio of equivalent stress multiplied by safety factor and of yield stress 
thnew(membrane#) - new skin thickness 
obliteration#) - mass after optimizing 
rmax - largest constant violation 

Following the optimization routine, the stresses are calculated again based on the new 
thickness. If the weight does not converge, then the routine returns to the optimization 
loop. 

A.4 Node and Element Numbering 

A description of the node and element numbering follows. The nodes are 
numbered on the upper surface then the lower surface. The first node is at the leading 
edge of thfroot chord. The nodes are then numbered along the chord to the trailing 
edge. Numbering continues in this manner for each chord from the root to the tip. 
Figure A1 illustrates this numbering scheme. Membranes are numbered in the same 
order as nodes. The upper skin is first followed by the lower skin. Rib webs are 
numbered next, then spar webs. Truss elements begin with spar caps on the upper, then 
lower surface. Vertical spacers are numbered next. Rib caps are then numbered, 
alternating from upper to lower surface for each rib section. (See figure Al) Buckling 
panels are numbered in the same order as membrane elements. 
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rib 



spar 


n# - indicates node number 
m# - indicates membrane number 
t# - indicates truss number 


Figure A 1. Node and Element Numbering 
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Point Loads: 

Point loads for fuel, engines, external stores, etc. Can be specified in a data File. 
The location of the point load must be specified within the geometry of the wing box. 
Shape functions are used to distribute the load to nodes of the cell which the point is 
located. Sum of forces and moments are conserved. 

Aerodynamic Loads: 

The aerodynamic loads which do not lie in the wing box are not included in the 
transformation of aerodynamic loads to structural points. These loads are placed at the 
nearest nodes along the leading or trailing edge. This method does not account for the 
change in moments. The change in moments could be dealt with by adding couple on 
the wing. 

Uniformly Distributed Load: 

CDOSS uses the gross weight estimate and the number of ‘g’s of the maneuver 
to calculate the distributed load. The gross weight multiplied my the number of ‘g’s is 
spread uniformly over the wing box area. These forces are then lumped onto the nodes 
based on the amount of area around each node. 



Natural Frequencies: 

To find natural frequencies of a wing, it is necessary to find the eigenvalues of 
the system. The general equation of motion is [-(D 2 M + K]{<|>} = {0} where co is the 
natural frequency, M is the mass matrix, and K is the stiffness matrix, and <J> is the 
vector of displacements. The eigenvalue solver solves equations of the form 
[A - A,I]{x} = {0} where I is the identity matrix. The equation of motion can be 
manipulated into this form. 

[-G) 2 M + K]{<t>} = {0} 

[- (O 2 4m 4m + Ar]{0} = {0} 

Va7 [- o 2 4 m 4 m + k]{<i>} = {o} 

[- co 2 4 m + 4mt'k]{(I)} = {0} 

let {<p} = 4 m~'{$} 

[- co 2 4 m + 4m^k\4m~ x {£} = {o} 

[- co 2 1 + 4m~'k4m~' ]{£} = {o} 

let X = to 2 

A = 4mFk4m ~ 1 

[A - XI]{£} = {0} 

The stiffness matrix is already formed for the stress analysis. The mass matrix is 
created by a simple lumped mass method. The weight of all structural elements are 
summed and lumped at nodes. No dynamic or external masses are included. The 
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mass matrix is a simple diagonal matrix. The inverse square root if this matrix is 
simply the reciprocal of the square root of the diagonal elements. The multiplication is 
carried out to obtain the A matrix. The eigenvalue solver finds the eigenvalues of A. 
To obtain natural frequencies in units of hertz, simply divide the square root of the 
eigenvalue by 2tz. 


(o{ hertz) = 77— 

2k 


The natural frequencies can be used for flutter analysis. 


Input and Output Files: 

Input Files (in units of feet, lbs, etc.): 

Geometry file - f 16 
gross weight 
number of components 
number of cross-sections 
number of points per cross-section 
list of points: x, y, and z coordinates 

Design data file - wingdata 
spar configuration 
number of spars 

number of ribs in section 1, section 2, etc. 
wing box shape 

percent chord for slats and flaps 
wingskin thickness at root and tip 
spar and rib web thickness 
spar and rib cap area 
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percent thickness for dummy webs 
vertical spacer area 
load choice 

aerodynamic transformation method or number of ‘g’s of maneuver 

filename of aerodynamic loads 

point loads, yes or no 

filename of point loads 

optimize, yes or no 

elastic modulus and density of bars 

minimum truss area 

maximum truss area 

elastic modulus, poisson’s ratio, and density of membranes 

minimum thickness of membranes 

maximum thickness of membranes 

step size for buckling optimization 

type of loading for membranes 

design safety factor 

material yield strength 

convergence criteria tolerance 

Aerodynamic Loads file - fl61oad 

• number of points 

• list of points: x, y, and z coordinates, vertical force 

Point Loads file - ptloads 
number of points 

list of points: x, y, and z coordinates, vertical force 

FEM files - wingtemp and opt.fl6 
filename 

number of nodes, truss elements, and membrane elements 
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list of nodes: node number x, y, z coordinates, x, y, z forces 
list of truss elements: element number, node 1, node 2, area 
list of membrane elements: element number, node 1 thru 3, thickness 

Displacement Output - out.fl6 

list of displacements: node number, x, y, z displacements 

list of truss element stresses: element number, stress 

list of membrane stresses: element number, principle stresses, shear stress 

Buckling Output - out.buckle 

list of panels: panel number, buckling constraint value 

list of truss elements: element number, buckling constraint value 

list of eigenvalues 
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